<!DOCTYPE HTML PUBLIC "-//W3C//DTD HTML 4.01 Transitional//EN" "http://www.w3.org/TR/1999/REC-html401-19991224/loose.dtd">
<html>
<head>

<title>Simulations for Statistical and Thermal Physics</title>

<link href="../default.css" type="text/css" rel="stylesheet">
</head>

<body>

<h3 style="text-align: center;">The Fermi-Pasta-Ulam problem</h3>

<p class="header_title">Introduction</p>

<p>The Fermi-Pasta-Ulam problem is named after the 
numerical experiments first performed by Enrico Fermi, John Pasta, 
and Stanislaw Ulam in 1954&#8211;1955 on the Los Alamos 
MANIAC computer, one of the first electronic computers. They wanted to understand how a crystal evolves toward thermal equilibrium by simulating a chain of particles coupled by springs.</p>

<p>&nbsp;&nbsp;&nbsp;&nbsp;Suppose that we start with a simple model of a crystal, a system of particles of mass m connected by harmonic springs &#8211; a model inspired by the Debye theory of a crystal. The energy of such a system can be expressed as:</p>
<p class="center">
<img src="eho.jpg" alt="" align="middle" >
</p><p>where u<sub>i</sub> is the displacement of particle i from its equilibrium position, p<sub>i</sub> is its momentum, and k is the usual spring constant. The two ends of the chain are assumed to be fixed and hence u<sub>0</sub> = u<sub>N</sub> = 0.</p>

<p>&nbsp;&nbsp;&nbsp;&nbsp;The usual way of understanding a system of masses connected by springs is to express the energy in terms of the normal modes. The latter are related to the displacement by</p>
<p class="center">
<img src="ft.jpg" alt="" align="middle" >
</p><p>
The energy of the system can be expressed in terms of the normal mode amplitudes as</p>
<p class="center">
<img src="enm.jpg" alt="" align="middle" >
</p><p>where</p>
<p class="center">
<img src="omegak.jpg" alt="" align="middle" >
</p>

<p>&nbsp;&nbsp;&nbsp;&nbsp;The equipartition theorem in statistical mechanics says that the average contribution of each quadratic contribution to the energy is kT/2, where T is the absolute temperature and k is Boltzmann's constant. However, the energy in each normal mode remains constant because the normal modes are independent of each other. Hence, if we solve Newton's equations of motion for a system of masses connected by linear springs, the energy in each mode will be the same as it was initially. That is, if we initially put all of the energy into a single mode, the energy will stay in this mode indefinitely, and the equipartition theorem will not be satisfied.</p>

<p>&nbsp;&nbsp;&nbsp;&nbsp;A system that is <i>ergodic</i> (or more accurately quasi-ergodic) is one for which all accessible microstates are equally probable over a long period of time. If the system were ergodic, the system would eventually have an equal 
distribution of energy between the various Fourier modes. Hence, we conclude that a system of particles connected by linear springs is non-ergodic.</p>

<p>&nbsp;&nbsp;&nbsp;&nbsp;What happens if we make the springs nonlinear? For example, suppose we take the potential energy to be</p>
<p class="center">
<img src="nonlinear.jpg" alt="" align="middle" >
</p>
<p>The last term leads to a coupling between the modes. Fermi, Pasta, and Ulam thought that this coupling would lead  to the energy that was  initially in a single mode eventually drifting into the other modes until the equipartition theorem is satisfied. However, when they did the simulation, they found to their great surprise that this energy was shared by only 
a few other modes; the remaining modes were hardly excited. Moreover, after a long time the initial state was almost completely recovered. This result, which is known as the FPU problem or the FPU paradox, shows that nonlinearity is not sufficient to guarantee the equipartition of energy. It turns out that ergodic behavior is only observed when the magnitude of the nonlinear term is more than a certain critical value.</p>

<center>
<applet
 code="org.opensourcephysics.davidson.applets.ApplicationApplet.class"
 archive="./stp.jar" codebase="../" align="top" height="40"
 hspace="0" vspace="0" width="150"> <param name="target"
 value="org.opensourcephysics.stp.fpu.FPUApp"> <param name="title"
 value="Applet"> <param name="singleapp" value="true">
</applet>
</center>

<p class="header_title">The Program</p>

<p>The program uses the Verlet algorithm to solve Newton's equations of motion numerically. The user can select the initial mode (k &#8805; 1) of the system. The energy in this mode (<font color = "red">red</font>) and the two neighboring modes k + 1 (<font color = "blue">blue</font>) and k - 1 (<font color = "green">green</font>) are plotted. The mean energy in these three modes is given when the simulation is stopped.</p>

<p>&nbsp;&nbsp;&nbsp;&nbsp;If the
system is ergodic, then all particles will see the same average
environment, and the time average of any physical quantity associated with individual particles
will be the same for each particle if
the time t is sufficiently long. The time average of this quantity f for particle i is defined as</p>
<p class="center">
<img src="fibar.jpg" alt="" align="middle" >
</p><p>The time 
average of f averaged
over all particles is</p>
<p class="center">
<img src="averageall.jpg" alt="" align="middle" >.
</p>

<p>&nbsp;&nbsp;&nbsp;&nbsp;We will choose the physical quantity of interest to be the potential energy of a
particle, e<sub>i</sub>, defined as</p>
<p class="center">
e<sub>i</sub> = &#931;<sub>i &#8800; j</sub>
u(r<sub>ij</sub>),
</p><p>where u(r<sub>ij</sub>) is the potential energy of interaction between particle i and j.</p>

<p>&nbsp;&nbsp;&nbsp;&nbsp;A measure of the difference between the mean value of e<sub>i</sub> compared to its average over all particles is the energy <i>metric</i>
&#937;(t), defined as</p>
<p class="center">
<img src="metric.jpg" alt="" align="middle" >.
</p><p>If the system is ergodic over the time
interval t, then it can be shown that
&#937;(t) decreases as
1/t. The program plots 1/&#937;(t) versus t.</p>

<p class="header_title">Problems</p>

<ol>

<li>First consider the behavior of a system of harmonic springs and set b = 0. Describe some of the normal modes for N = 4. Does the energy that is initially in a given mode change with time? What is the behavior of the metric?</li>

<li>Use the default parameters, but choose b = 0.25. What is the behavior of the metric for t &#60; 2 &#215; 10<sup>5</sup>? What is the behavior of the energy in the three modes that are plotted? What is the behavior of the energies and metric for longer times? (The unit of time is such that k = m = 1.) The simulation can be speeded up by choosing the number of steps per display to be at least 1000 and closing the window that displays the displacements of the particles.</li>

<li>Increase the value of b and determine the approximate value of b for which the energy in the various modes becomes approximately equal and the inverse metric increases approximately linearly.</li>
</ol>

<p class="header_title">References</p>

<ul>

<li>David K. Campbell,
Phillip Rosenau, and
George M. Zaslavsky, "The Fermi-Pasta-Ulam problem&#8211;The first fifty years,"" Chaos <b>15</b>, 015101 (2005).</li>

<li>Thierry Dauxois, Michel Peyrard, and Stefano Ruffo, "The Fermi-Pasta-Ulam 'numerical 
experiment': History and pedagogical 
perspectives," Eur. J. Phys. <b>26</b>, S3&#8211;S11 (2005).</li>

<li>Thierry Dauxois and Michel Peyrard, <i>Physics of Solitons</i>, Cambridge University Press (2006), Chapter 8.  Parts of this graduate level text are acessible to undergraduates.</li>

<li>T. P. Weissert, <i>The Genesis of Simulation in Dynamics: Pursuing the Fermi-Pasta-Ulam Problem,</i> 
Springer (1997).</li> 

</ul>

<p class = "small">Updated 2 May 2007.</p>

</body>
</html>


